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Abstract 

This paper considers new approach to control a stationary inhomogeneous planar 
front solution of a nonlinear parabolic two-dimensional distributed (reaction-diffusion) 
system, by using a gain point-sensor control with actuators that have the simplest 
possible spatial dependence. The method is based on multivariable root-locus tech- 
nique for the finite-dimensional approximation of the original PDE model and use the 
concepts of finite and infinite zeros of linear multidimensional system. 
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1 Introduction 

Nonlinear parabolic partial differential equations (PDEs), which typically describe reaction- 
diffusion systems, may admit spatially-dependent solutions like stationary fronts as well as 
spatiotemporal patterns. The latter can often be described as composed of slow-moving 
fronts, separated by domains of moderate changes. This article is part of a research aimed to 
develop control theory for reaction-diffusion and reaction-convection-diffusion systems for 
which a certain patterned state is advantageous. Propagating fronts and patterned states 
may emerge in several technologies including catalytic reactors [1], distillation processes 
[2], flame propagation and crystal growth [3] (see also [4] for references) as well as in 
physiological systems like the heart [5]. Our interest lies in catalytic reactors, in which 
stationary or moving fronts and spatiotemporal patterns have been observed and simulated 
in various systems like flow through a catalyst [6], fixed-bed reactors [7, 8], reactors with 
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flow reversal [9] and loop reactors [10]. The instabilities emerge due to the thermal effects 
in exothermic reactions, due to self-inhibition by a reactant and due to slow reversible 
modifications of the surface. The construction of a controller that enables to stabilize 
some inhomogeneous solutions in one-dimensional (1-D) reaction-diffusion and reaction- 
convection-diffusion systems is currently a subject of intensive investigation [11, 12, 13, 14, 
15, 16, 17]. Yet, most catalytic reactors, as well as physiological systems like the heart, 
exhibit a behaviour that can be described by two- or even three-dimensional models. 

In the present work we are interested in stabilizing stationary (planar) fronts in a rectan- 
gular two-dimensional (2-D) domain in which a chemical reaction-diffusion process occurs. 
Previous studies of 1-D systems demonstrated that the simplest approach is by applying 
the point-sensor control, in which a single space-independent actuator responds to a sensor 
that is located at the front position [16] . In the 2-D problem the point-sensor control with a 
single space-independent actuator cannot stabilize the front in a wide system. In that case 
it is necessary to use control with several point sensors and actuators (space-independent 
and space-dependent). Here we apply the extension of the point-sensor control to multi- 
variable case by explore the root-locus technique for multivariable system [21]. The main 
advantage of this control is its insensitivity to parameter uncertainties (robustness) and also 
small number of evaluated gain coefficients. We apply the root-locus method to determine 
the minimal number of actuators , their spatial form and the location of the correspond- 
ing sensors that will assure the linear stability of the planar front. The method uses the 
concepts of finite [20, 19] and infinite zeros [21] of a linear multidimensional system; these 
notions are briefly explained below. The proposed approach is best suited for a systematic 
computer-aided search of the regulator form. 

2 Statement 

Consider the reaction-diffusion problem, in the (z, r) rectangular domain of length L and 
width R , which is described by a pair of coupled nonlinear parabolic PDEs 

y t -y zz -y r r = P(y,0) + \ 9 t = eQ(y,6) (1) 
subject to no-flux boundary conditions : 

y 2 (0,r) = 0, y z (L,r)=0, y r (z, 0) = 0, y r (z,R) = (2) 

where the variable y = y(z,r,t) typically represents the activator, 9 = 0(z,r,t) is the slow 
variable (localized inhibitor), A is a control variable that is introduced additively; e < 1 is 



the ratio of time scales. 

We use the polynomial source functions 



P( y ,9) = -y ;i + y + 9, Q(y,0) = ->yy - 



(3) 



since they adequately simulate the phenomenon of multiple steady states (bistable kinetics) 
due to thermal and autocatalytic effects which, in many chemical systems, induce the 
instabilities [22] and since several analytical results are available [23] . The instability stems 
from the anticlinal arrangement of y and 9 [4]. In Eqn.(3) 7 > is a constant parameter. 
System (1-3) has often been used for studying pattern-formation in chemically reacting 
systems [1]; similar systems were employed to describe pattern formation and control in 
cardiac systems. 

In this work we address the problem of stabilizing the planar front, in a 2-D reaction- 
diffusion system (1-3), at the position which in an open- loop system is stationary but 
unstable (i.e., at the middle, z = L/2 of the domain). Obviously, narrow systems behave 
like a one-dimensional system (Eqn. (1-3), y rr = 0) and admit a stationary front solution 
that is unstable and typically oscillates or travels out of the system (see details in [17]). In 
a sufficiently wide 2-D system, the front-line may undergo symmetry breaking so that in 
part of it the upper state expands while in other parts the lower state propagates. In that 
case we suggest to use control with several actuators located along the front at the points 
(z,r) = (L/2, r d ), d = 1, . . . , 77 and apply a general feedback control law of the form 



where k < is a scalar gain coefficient, y(L/2, r d , t) — y* d are deviations of the sensors from 
the set points; y* d = y*(L/2,r d ); ip d (z,r) are some space-dependent functions that may 
imitate the eigenfunctions. We will seek control (4) with the simplest space-independent or 
space-dependent actuator functions ip d (z,r) and minimal number of sensors. 

3 Linear analysis 

For design of control (4) we use a linearized truncated version of PDEs (1-3). Linearization 
of (1-3) for A = around the steady state solution y s = y s (z,r), 9 S = 9 s (z,r) yields 
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\ = kY,[y(L/2,r d ,t)-y d }M^r) 



(4) 



d=l 



yt - Vzz - y^ = 



(-3f s + l)y + 9 + \, 



9 t = -e>yy - e9 



(5) 



which is lumped by the Galerkin method: We expand the deviations y = y(z,r,t), 9 = 
6(z, r, t) and y(z*, r d , t) = y(z*, r d , t) - y s (z*, r d ) as 



y( z , t) = Yl a e{t)4>e{z, r), 9(z, t) = b e (t)<j>e(z, r), y(z* , r d , t) = J2 a e (t)(f) e (z* , r d ) (6) 

e e e 

where z* = L/2 and the functions (f> e (z,r) are the eigenfunctions of the problem 

<t>zz{z, r) + (f) rr (z, r) = -\</>(z, r), 2 (O, r) = <j> z (L, r) = <f> r (z, 0) = <j> r (z, R) = (7) 
with the corresponding eigenvalues 

\ em = A{i-l) 2 /L 2 + U-^?/R 2 } (8) 

and the eigenfunctions 

, p (i — l)nz (j — l)nz 

(Pe(ij) = -n=^ cos F cos 5 ' e = e w) = 2 ' • • • ( 9 ) 

y LK L, Li 

(p = 1 when % — j — 1, p — 2 when i,j>l and p = a/2 when i = 1, j > 1 or j = l,i > 1, 
see [4] for derivation ). Substituting (6) into (5), (4) with set points y d = y s (L/2,r d ) and 
integrating with a weight eigenfunctions <f> e (z,r) results in the spectral representation of 
closed- loop linearized system (5), (4) 

■n r-L rR 

a e = -\ e a e + Je m a m + K + k / / ^d(z,r)J2 a f ( Pf(ki)( L / 2 , r d)<Pe(ki)(z,r)dzdr (10) 

m d= \ Jo Jo f 

b e = e(7a e + b e ) (11) 

f L f R 

Jem = J o J Q (-tys + l)<p m (ij)4> e (ki)dzdr e, m = 1, . . . 



where 



Denoting 



h df = <f) f {L/2,r d ) 

fL rR 

(3ed= I I ip d {z,r)<f) eikl) (z,r)dzdr 



o Jo 



we can rewrite (10) as follows 



h e = -\ e a e +J2 J e™ a rn + K + kj^ PedJ2 ° f h <V ' e > / = 2 > • • • 
m d=l f 



(12) 

(13) 
(14) 

(15) 



This system (Eqn. 15, 11) may be presented in the usual vector-matrix form as the linear 
infinite-dimensional dynamical system with 77 - dimensional input v and output w vectors 



(16) 
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w = Ha 

closed by the finite-dimensional output feedback 



(17) 



v = kl v w (18) 

where a(t) = [a e ],b(t) = [b e ] are the infinite dimensional vectors (e = 1,2,...); v and w 
are finite-dimensional rj vectors, the matrix f3 = \J3 e( j\ has 77 infinite-dimensional columns 
(e = 1,2, ...;d = 1, . . . , 77) and the matrix H = [hjf] has 77 infinite-dimensional rows 
(d — 1, . . . , 77; / = 1,2,...); I v is unity 77 x 77 matrix; A = diag(\i, X 2 , . . .), I — diag(l, 1, . . .) 
, J = [J em ], e, m = 1, 2, . . . are infinite-dimension matrices and k is scalar gain coefficient. 

Let us study the correlation between the form of original control (4) and the input and 
output structures of Eqns. (16-17): If all r d , d — 1, 2, . . . , 77 in (4) have different values, then 
the parameter 77 assigns the dimension of input and output vectors of Eqns. (16), (17) and 
the sensor positions stipulate the structure of the output matrix H (see Eqn.13). Besides, 
the form of the actuator functions ipd{z,r) influences the matrix (3 (see Eqn 14). Hence, in 
view of the original statement the problem can be stated as follows: 

Problem 1. For the linearized infinite-dimensional ODEs system (16-17) it is necessary 
to find the matrices f3 and H with the minimal number of columns and rows, respectively, 
such that the finite-dimensional output feedback control (18) stabilizes the closed-loop 
system. 

The obvious way of designing of a finite-dimensional control (18) for an infinite-dimensional 
system (16-17) is to use a truncated (finite-dimensional) approximation of the PDEs. To 
truncate the system we capitalize on the dissipative nature of the parabolic PDEs: the trun- 
cation order N is estimated by calculating the leading eigenvalues of the dynamics matrix 
(16) (see [4] for details). As follows from [24], for a sufficiently large truncated order N, with 
actuator functions ipd(z, r) that coincide with eigenfunctions (f> e (z, r), the finite-dimensional 
control guarantees the stability of infinite-dimensional system. Below we imply that system 
(16-17) is the finite-dimensional analog of the original (16-17) with truncated order N. 

4 Root-locus control design 

Design of control (4) implies the determination of the minimal number of actuators to be 
employed, their spatial form and the location of the corresponding sensors that will assure 
the linear stability of the planar front. We can manipulate the number and r-coordinate 
positions (ri, r 2 , . . . , r v ) of the sensors and the form of the actuator functions ipd(z, r )- The 
former affects the matrix H (Eqn. 13) and the latter influences the matrix (3 (Eqn. 14). 



For simplification of the search for a matrix (3 structure we assign the eigenfunctions (9) 
as actuator functions, i.e. ipd( z , r ) ~ <j>e( z i r),e — 1,2, . . .. Thus from the relation between 
the form of the actuator functions ipd{z, r) and the structure of the matrix [3(see Eqn.14) it 
follows that every d-th column of the matrix (3 will contain only a single non-zero element 
f3 e d ■ Moreover, such a structure of f3 satisfies the above-mentioned Balas's restrictions [24] 
that ensure the validity of control of an infinite-dimensional system by a finite- dimensional 
controller. 

Therefore, we propose the following steps for design of control (4): (i) Assign sensor 
positions along the front (z,r) = (L/2,r^), d = 1, ... 77 and calculate the matrix H. 
(ii) Seek the matrix (3 and gain k that assures linear stability of the closed-loop system 
(Eqns. 16-18). (iii) Finally, find the eigenfunctions <f> e (z,r) that corresponds to the matrix 
(3 obtained. 

For fulfilling step (ii) we use a root-locus technique [21] which is based on an analysis 
of the finite zeros and infinite zeros of the open-loop system (16-17). This method uses 
the following known property of closed-loop linear system with feedback (18) (see [25]): as 
the feedback gain increases towards infinity a part of the closed-loop eigenvalues remain 
finite and approaches the values which are referred to as finite system zeros [19], [20] (see 
also Appendix for definitions) while the remainder are located at the points at infinity and 
are known as infinite zeros [21]. Therefore, we propose to seek a suitable matrix f3 by 
the repeatedly calculating the finite zeros of open-loop ODEs (16-17), with different input 
matrices, and finding the one that ensures that the leading finite system zeros are negative. 
Then we assign a sufficienly large negative gain coefficient k. 1 

Finally, Problem 1 may be reformulated as follows: 

Problem 2. For the linearized truncated system (16-17) with assigned output matrix 
H it is necessary to find the matrix (3 that has a single nonzero element in every column 
such that the above-mentioned system has leading finite system zeros in the left -half of the 
complex plane ('negative' zeros) 2 . 

Let us note that such input and output matrices H and f3 result a minimum phase 
control system because this system contains finite zeros in the left-half of the complex 
plane. 

The following assertion is needed to avoid cases when Problem 2 has no solution . 

1 Here we must note that large perturbations of the front will cause the control variable to exceed the 
bistability domain of P = . So the gain value coefficient k must be restricted below the value (k ) that 
ensures stability of the closed-loop system. 

2 Similarly we denote zeros in the right-half of the complex plane by 'positive' ones. 



Assertion 1. For assiged matrix H the problem has no solution for any matrix (3 if 
and only if 'positive' finite system zeros of (16), (17) are output-decoupling zeros [26]. 

The proof follows from definition of decoupling zeros (see Appendix). 

Remark 1. If the shape of actuator distribution functions ipd{z,r),d = 1,2,... are 
preassigned from technical constraints (i.e. the matrix f3 is given) then we need to find 
r] sensor locations (i.e. the matrix H ) which provide the 'negative' leading finite system 
zeros of system (16-17). The solvability of this problem is formulated as follows. 

Assertion 2. For assiged matrix (3 the problem has no solution for any matrix H if 
and only if 'positive' finite system zeros of (16), (17) are input-decoupling zeros [26] . 

Remark 2. It is necessary to choose the sensor locations in r-direction so that the 
infinitely increasing eigenvalues of the closed-loop high gain system (16-18) (infinite zeros) 
tend to infinity along asymptotics with a negative real angle. This condition is guaranteed 
if det(H(3 7^)0 and all eigenvalues of the r\ x r\ matrix Hf3 are positive [21]. If H(3 has 
several negative eigenvalues then we need to introduce a nonsingular precompensator M to 
(18) such that the new control 



ensures above property for MH(3 . Such operation does not changes the finite system 
zeros of (16), (17) because they are invariant to any nonsingular transformation of output 



Therefore, the general strategy of the method is as follows: At first we need to check that 
assigned sensor positions (or matrix H) ensure that output-decoupling zeros of a pair (A, H) 
are 'negative' (see Assertion 1). Then we seek a suitable matrix (3 by repeating calculations 
of the finite system zeros of open-loop ODEs (16), (17) with different matrices (3 and 
finding one that ensures that the leading finite system zeros are negative and det(H[3) ^ 0. 
If necessary, we find the precompensator M which rearranges the sensor position in re- 
direction so that the infinite zeros tend to infinity along asymptotes with a negative real 
axis angle. 

To demonstrate this procedure we apply it for design of control (4) for PDEs (1-3) of 
length L and various widths. 

We start by analyzing the effectiveness of the simplest control law, a single space- 
independent actuator (Eqn. 4, r\ = 1, ipi(z,r) = 1) 



v = kI v Mw 



(19) 



[19], [20]. 



\(t)=ky((L/2, ri ,t)-y* a ) 



(20) 



where y* = y s (L/2,ri) coincides with the steady-state value of the problem . The spectral 
representation of the closed-loop system (1-3), (20) is a single-input, single-output Eqns.(16- 



17) with a column vector (3 = [/3n, 0, . . .] with fin = 1/yLR and row output vector H = 
h = [hi, h 2 , . . .]. Thus here the shape of the actuator is assigned. Then if input-decoupling 
zeros of a pair (A, (3) are 'negative' we need to find a sensor position. At first consider a 
sensor situated at the domain center (z* = L/2,r\ = R/2). The analysis of leading system 
zeros of the related linearized truncated system obtained shows that they possess negative 
real parts for < R < R cr for some value R cr and positive real parts for R > R cr . So, 
control with one space-independent actuator is effective only for systems of width R < R cr 
(narrow systems) that have 'negative' zeros. Changing the sensor position in the r -direction 
alters the vector h (Eq.13 ) and as a consequence may move the finite system zeros (for 
details see [27]) and R cr . However, the position at the domain center assures the maximal 
value of R cr . 

For wider systems (R > R cr ) it is necessary to use space-dependent actuators. We try 
to apply control in the form 

A = k^2[y(L/2,r d ,t)-y:]Mz,r) (21) 
d=i 

with one space-independent actuator (^i(z, r) = 1) and another space-dependent one, 
il) 2 (z,r) ~ (f) e (z,r) where <p e (z,r) is eigenfunction (9) chosen from the series eigenfunc- 
tions ordered in an increasing order of the appropriate eigenvalues. Introducing two sensors 
at positions (L/2,r 1 ) and (L/2,r 2 ) we calculate the 2 x N matrix H by Eqn. (13). If 
output-decoupling zeros of a pair (A, H) are 'negative' then we need to evaluate the fi- 
nite system zeros of systems (16-17) with above output matrix H and different N x 2 
matrices (3 = [f3 e d], e = l,...,N,d = 1,2 with assigned elements of the first column 
(Pu = 1/ V LR, P e i = 0, e = 2, • • • , N ) and undetermined elements of the second column. 

Thus it is necessary to find a single nonzero element from the second column of (3 which 
ensures that the leading finite zeros of system (16-17) are 'negative' ones. If such (3 does not 
exist we need to change the matrix H by shifting the positions of sensors in r-direction and 
begin the search of (3 once more. If an appropriate matrix H does not exist then control 
(21) is not effective and it is necessary to increase the number of sensors (and actuators) 
to three and so on. 

5 Application 

Let us apply the above method for stabilization of planar front solution of system (1-3) 
(with L = 20, 7 = 0.45, e = 0.1 ) of various widths. The leading eigenvalues of the 
truncated version (Eqn. 16, iV = 23) show two real unstable eigenvalues ( 0.35 and 0) for 



all i?'s and two complex eigenvalues with real parts that becomes positive for R > 5.6( see 
Fig. 4a in [4]). At first we try to apply control (20) with one space-independent actuator 
and the sensor at (z* = 10, r\ = R/2). The analysis of zeros of open-loop system (16-17) 
discovers two leading (and complex) finite zeros with negative real parts for < R < 5.6 
and positive real parts for R > 5.6 (see Fig. 4b in [4]). Hence, system width R cr = 5.5 is a 
critical one for our ability to do control with a space-independent actuator. 

For wider systems (R > 5.5) we apply the two actuator control (21). Assigning two 
sensors in some positions (z* = 10, ri), (z* = 10, r 2 ) (n ^ r 2 ) and calculating the finite 
zeros of relevant two input, two-output systems (16-17) with different input matrices we 
find that when (3 24 7^ the leading finite zeros have negative real parts (—0.1024). This 
correspondents to the fourth (e = 4) eigenfunction 4 ~ cos(nr/R) in series of the ordered 
eigenfunctions : The first six ordered eigenfunctions are 0i = 1, <p 2 ~ cos(nz/ L), 03 ~ 
cos(27iz/L), 04 ~ cos{nr / R), 5 ~ cos{?>Trz/ L), 6 ~ cos(2iiz / L)cos(nr / R) . Consequently 
control (21) becomes 

X(t) = k{[y(10, ri ,t) -y s (10, ri )\ + [y(10,r 2 ,t) - y s (10,r 2 )]cos(7rr/R} (22) 

Then, it is necessary to choose r±, r 2 in (22) so that Remark 2 is satisfied. For our case 
we need to use the sensor locations with r\ > r 2 which is equivalent to introducing the 2x2 

1 

1 

We verify the proposed methodology by simulating system (1-3) using one space- inde- 
pendent actuator (20) for narrow (R < R cr ) and two-actuator control (22) for wider system 
(see [4] for plots). 



precompensator M 



6 Conclusion remarks 

The stabilization of planar stationary fronts in a two-dimensional rectangular domain, in 
which a diffusion-reaction systems occurs, is studied using a two-variable PDEs model for 
which some analytical results are available. We consider the simplest control strategy based 
on sensors placed at the designed front line position and measure deviations from a local 
state, and actuators that are spatially-uniform or space dependent. We present a systematic 
control design that determines the number of required sensors and actuators, their position 
and their form. The control design is corroborated by linear analysis of a lumped truncated 
model and concepts of finite and infinite zeros of linear multidimensional systems. The 
method is best suited for a systematic computer-aided search of the regulator form. 



APPENDIX 



Consider a general linear multivariable finite-dimensional dynamic system described by 
the set of state-space equations 

d(f) = Aa(t) + Bv(t) w(t) = Ha(t) 

with a , v and w are the n , r and / dimensional state, input and output vectors and A, B 
and H are constant matrices of appropriate dimensions. 

Definition 1. The finite zeros (system zeros) of the above system are determined as 
the set of complex s; for which the rank of the system matrix 



P(s) 



Sil-A -B 
H O 



is reduced. 

Definition 2. The input-decoupling zeros are defined as the set of the complex variable 
Si at which the row rank of the matrix [sj — A, —B] is reduced. 

Definition 3. The output-decoupling zeros are defined as the set of the complex variable 



Si at which the row rank of the matrix 



is reduced. 



Si I - A T , H T 

Evaluation of input-decoupling zeros. The set of input-decoupling zeros may be 
calculated from Definition 2. To avoid the operations with complex numbers we can apply 
an alternative method which uses the property: input-decoupling zeros coincide with the 
uncontrollable eigenvalues of the matrix A. Since the latter eigenvalues are invariant under 
a proportional state feedback: v = Ka then they may be calculated as those eigenvalues 
of the closed-loop matrices A + BKj which are invariant with respect to any gain matrices 
Kj with finite elements. Let us analyze the leading eigenvalues of the iV-truncated system 
(16) acted by the state control v = Kja with a random iV row vector Kj. It is evident that 
leading eigenvalues, which are invariant with respect to this control, are input-decoupling 
zeros of (16), (17). 
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